*
* $Id$
*
* $Log: photon.F,v $
* Revision 1.1.1.1  2002/07/24 15:56:27  rdm
* initial import into CVS
*
* Revision 1.1.1.1  2002/06/16 15:18:43  hristov
* Separate distribution  of Geant3
*
* Revision 1.1.1.1  1999/05/18 15:55:21  fca
* AliRoot sources
*
* Revision 1.1.1.1  1995/10/24 10:21:58  cernlib
* Geant
*
*
#include "geant321/pilot.h"
*CMZ :  3.21/04 23/02/95  14.46.01  by  S.Giani
*-- Author :
      SUBROUTINE PHOTON(D,LD,IDICTS,LDICT,NTX,NTS,IGCBS,LGCB,
     +      AWR,IGCBS2,LGCB2,LR,IGAMS,LGAM,QI,ID,IIN,LRI,SIGN)
C       THIS ROUTINE CONTROLS THE GENERATION AND STORAGE OF ALL
C       PHOTONS PRODUCED BY THE NEUTRON INTERACTIONS.  WHERE DATA
C       PERMITS, THE PHOTON PRODUCED IS DIRECTLY COUPLED TO THE
C       NEUTRON REACTION OCCURING.
#include "geant321/minput.inc"
#include "geant321/mconst.inc"
#include "geant321/mnutrn.inc"
#include "geant321/mapoll.inc"
#include "geant321/mcross.inc"
#include "geant321/mpstor.inc"
#include "geant321/mmicab.inc"
      DIMENSION IDICTS(NNR,NNUC),LDICT(NNR,NNUC),NTX(*),NTS(*),
     +  IGCBS(NGR,NNUC),LGCB(NGR,NNUC),AWR(*),IGCBS2(NGR,NNUC),
     +  LGCB2(NGR,NNUC),LR(NQ,NNUC),IGAMS(*),LGAM(*),D(*),LD(*)
      SAVE
C flag to mark call to SECEGY = 1 or PARTXS = 2 for EP   CZ 13/8/92
      IEP = 0
C       INITIALIZE THE PHOTON ENERGY TO ZERO IN CASE NO PHOTON IS
C       CHOSEN (THIS IS NECESSARY BECAUSE OF ENDF INCONSISTENCY)
      EG=0.0
C       INITIALIZE THE PARAMETERS USED IN THE SELECTION PROCESS
      MT=0
      IMT=0
      NUMBG=0
      XSIG2=0.0
      XSIG=0.0
      SIGMT3=0.0
      SIGP=0.0
      AWRI=AWR(IIN)
      NNTX=NTX(IIN)
      NNTS=NTS(IIN)
      L=2*NNTX+2*NNTS
C       NO PHOTON DATA PRESENT (IF L=0)
      IF(L.EQ.0)GO TO 360
      LX=2*NNTX
      LS=LX+1
C       DETERMINE THE NEUTRON REACTION MT NUMBER
      IF(ID.EQ.8)MT=16
      IF(ID.EQ.9)MT=17
      IF(ID.EQ.10)MT=18
      IF(ID.EQ.11)MT=22
      IF(ID.EQ.12)MT=24
      IF(ID.EQ.13)MT=28
      IF((ID.GE.14).AND.(ID.LE.54))MT=51
      IF(ID.EQ.55)MT=102
      IF(ID.EQ.56)MT=103
      IF(ID.EQ.57)MT=104
      IF(ID.EQ.58)MT=105
      IF(ID.EQ.59)MT=106
      IF(ID.EQ.60)MT=107
      IF(ID.EQ.61)MT=108
      IF(ID.EQ.62)MT=109
      IF(ID.EQ.63)MT=111
      IF(ID.EQ.64)MT=112
      IF(ID.EQ.65)MT=113
      IF(ID.EQ.66)MT=114
C       DETERMINE WHICH DISCRETE INELASTIC SCATTERING LEVEL OCCURRED
      IF(MT.NE.51)GO TO 130
      IMT=ID-14
      MT=MT+IMT
C       RESET THE MT NUMBER IF AN LR-FLAG IS INVOLKED
      IF(LRI.EQ.22)MT=22
      IF(LRI.EQ.23)MT=23
      IF(LRI.EQ.28)MT=28
C       CHECK PHOTON PRODUCTION DICTIONARY TO SEE IF THERE IS PHOTON
C       DATA CORRESPONDING TO THE NEUTRON MT REACTION THAT OCCURRED
      DO 10 IX=1,NNTX
         MTG=LGCB(2*IX-1,IIN)
         IF(MTG.EQ.MT)GO TO 30
   10 CONTINUE
   20 IF(LRI.EQ.22)GO TO 190
      IF(LRI.EQ.23)GO TO 190
      IF(LRI.EQ.28)GO TO 190
      GO TO 70
C       PHOTON DATA FOUND CORRESPONDING TO NEUTRON MT REACTION
   30 L1=LGCB2(2*IX,IIN)
      IF(L1.EQ.0)GO TO 370
      LS1=IGCBS2(2*IX,IIN)+LMOX4
      LEN=L1/2
      CALL TBSPLT(D(LS1),EOLD,LEN,SIGP)
      IF(SIGP.EQ.0.0)GO TO 190
      LS2=IGCBS(2*IX,IIN)+LMOX2
C       DETERMINE EXIT PHOTON ENERGY (EP)
      CALL PARTXS(D(LS2),D(LS2),EOLD,SIGP,EP)
      IEP = 2
      IF(EP.GT.0.0)GO TO 60
C       DISCRETE PHOTON ENERGY WAS NOT SELECTED (EP=0.0)
C       CHECK SECONDARY PHOTON DISTRIBUTION (FILE 15) FOR EP
      DO 40 IS=1,NNTS
         MTGS=LGCB(LX+2*IS-1,IIN)
         IF(MTGS.EQ.MT)GO TO 50
   40 CONTINUE
C no file 15 found and EP=0 in PARTXS -> try MT=4 etc
      GO TO 20
   50 L1=LGCB(LX+2*IS,IIN)
      IF(L1.EQ.0)GO TO 380
      LS3=IGCBS(LX+2*IS,IIN)+LMOX2
C       DETERMINE EXIT PHOTON ENERGY (EP)
      CALL SECEGY(EP,D(LS3),EOLD,D(LS3))
      IEP = 1
C       DETERMINE THE PHOTON MULTIPLICITY (YP)
C       RECALCULATE THE DENOMINATOR USED IN CALCULATING THE
C       PHOTON MULTIPLICITY TO ACCOUNT FOR THE LR-FLAGS
   60 IF(LRI.EQ.22)CALL LRNORM(D,D,IDICTS,LDICT,LR,EOLD,MT,IIN,SIGN)
      IF(LRI.EQ.23)CALL LRNORM(D,D,IDICTS,LDICT,LR,EOLD,MT,IIN,SIGN)
      IF(LRI.EQ.28)CALL LRNORM(D,D,IDICTS,LDICT,LR,EOLD,MT,IIN,SIGN)
      YP=SIGP/SIGN
      GO TO 330
C       THE DISCRETE INELASTIC LEVEL PHOTON DATA WAS NOT FOUND
C       CHECK THE PHOTON PRODUCTION DICTIONARY TO SEE IF THERE IS
C       PHOTON DATA CORRESPONDING TO MT=4
   70 DO 80 IX=1,NNTX
         MTG=LGCB(2*IX-1,IIN)
         IF(MTG.EQ.4)GO TO 90
   80 CONTINUE
      GO TO 190
C       PHOTON DATA FOUND CORRESPONDING TO MT=4
   90 L1=LGCB2(2*IX,IIN)
      IF(L1.EQ.0)GO TO 370
      LS1=IGCBS2(2*IX,IIN)+LMOX4
      LEN=L1/2
      CALL TBSPLT(D(LS1),EOLD,LEN,SIGP)
      IF(SIGP.EQ.0.0)GO TO 190
      LS2=IGCBS(2*IX,IIN)+LMOX2
C       DETERMINE EXIT PHOTON ENERGY (EP)
      CALL PARTXS(D(LS2),D(LS2),EOLD,SIGP,EP)
      IEP = 2
      IF(EP.GT.0.0)GO TO 120
C       DISCRETE PHOTON ENERGY WAS NOT SELECTED (EP=0.0)
C       CHECK SECONDARY PHOTON DISTRIBUTION (FILE 15) FOR EP
      DO 100 IS=1,NNTS
         MTGS=LGCB(LX+2*IS-1,IIN)
         IF(MTGS.EQ.4)GO TO 110
  100 CONTINUE
      GO TO 380
  110 L1=LGCB(LX+2*IS,IIN)
      IF(L1.EQ.0)GO TO 380
      LS3=IGCBS(LX+2*IS,IIN)+LMOX2
C       DETERMINE EXIT PHOTON ENERGY (EP)
      CALL SECEGY(EP,D(LS3),EOLD,D(LS3))
      IEP = 1
C       DETERMINE THE PHOTON MULTIPLICITY (YP)
C       RECALCULATE THE DENOMINATOR USED IN CALCULATING THE
C       PHOTON MULTIPLICITY TO ACCOUNT FOR THE LR-FLAGS
  120 MT=4
      CALL LRNORM(D,D,IDICTS,LDICT,LR,EOLD,MT,IIN,SIGNIS)
      YP=SIGP/SIGNIS
      GO TO 330
C       CHECK PHOTON PRODUCTION DICTIONARY TO SEE IF THERE IS PHOTON
C       DATA CORRESPONDING TO THE NEUTRON MT REACTION THAT OCCURRED
  130 DO 140 IX=1,NNTX
         MTG=LGCB(2*IX-1,IIN)
         IF(MTG.EQ.MT)GO TO 150
  140 CONTINUE
      GO TO 190
C       PHOTON DATA FOUND CORRESPONDING TO NEUTRON MT REACTION
  150 L1=LGCB2(2*IX,IIN)
      IF(L1.EQ.0)GO TO 370
      LS1=IGCBS2(2*IX,IIN)+LMOX4
      LEN=L1/2
      CALL TBSPLT(D(LS1),EOLD,LEN,SIGP)
      IF(SIGP.EQ.0.0)GO TO 190
      LS2=IGCBS(2*IX,IIN)+LMOX2
C       DETERMINE EXIT PHOTON ENERGY (EP)
      CALL PARTXS(D(LS2),D(LS2),EOLD,SIGP,EP)
      IEP = 2
      IF(EP.GT.0.0)GO TO 180
C       DISCRETE PHOTON ENERGY WAS NOT SELECTED (EP=0.0)
C       CHECK SECONDARY PHOTON DISTRIBUTION (FILE 15) FOR EP
      DO 160 IS=1,NNTS
         MTGS=LGCB(LX+2*IS-1,IIN)
         IF(MTGS.EQ.MT)GO TO 170
  160 CONTINUE
      GO TO 380
  170 L1=LGCB(LX+2*IS,IIN)
      IF(L1.EQ.0)GO TO 380
      LS3=IGCBS(LX+2*IS,IIN)+LMOX2
C       DETERMINE EXIT PHOTON ENERGY (EP)
      CALL SECEGY(EP,D(LS3),EOLD,D(LS3))
      IEP = 1
C       DETERMINE THE PHOTON MULTIPLICITY (YP)
  180 YP=SIGP/SIGN
      GO TO 330
C       NO PHOTON DATA WAS FOUND FOR THE PARTICULAR NEUTRON MT
C       REACTION OR FOR NEUTRON MT=4, THEREFORE CHECK THE PHOTON
C       PRODUCTION DICTIONARY TO SEE IF THERE IS PHOTON DATA
C       CORRESPONDING TO MT=3 (THE CATCH-ALL MT)
  190 DO 200 IX=1,NNTX
         MTG=LGCB(2*IX-1,IIN)
         IF(MTG.EQ.3)GO TO 210
  200 CONTINUE
      GO TO 360
C       PHOTON DATA FOUND CORRESPONDING TO MT=3
  210 L1=LGCB2(2*IX,IIN)
      IF(L1.EQ.0)GO TO 370
      LS1=IGCBS2(2*IX,IIN)+LMOX4
      LEN=L1/2
      CALL TBSPLT(D(LS1),EOLD,LEN,SIGP)
      IF(SIGP.EQ.0.0)GO TO 360
      LS2=IGCBS(2*IX,IIN)+LMOX2
C       DETERMINE EXIT PHOTON ENERGY (EP)
      CALL PARTXS(D(LS2),D(LS2),EOLD,SIGP,EP)
      IEP = 2
      IF(EP.GT.0.0)GO TO 240
C       DISCRETE PHOTON ENERGY WAS NOT SELECTED (EP=0.0)
C       CHECK SECONDARY PHOTON DISTRIBUTION (FILE 15) FOR EP
      DO 220 IS=1,NNTS
         MTGS=LGCB(LX+2*IS-1,IIN)
         IF(MTGS.EQ.3)GO TO 230
  220 CONTINUE
      GO TO 380
  230 L1=LGCB(LX+2*IS,IIN)
      IF(L1.EQ.0)GO TO 380
      LS3=IGCBS(LX+2*IS,IIN)+LMOX2
C       DETERMINE EXIT PHOTON ENERGY (EP)
      CALL SECEGY(EP,D(LS3),EOLD,D(LS3))
      IEP = 1
C       THE PHOTON WAS SELECTED FROM PHOTON DATA FOR MT=3
C       TO OBTAIN THE CORRECT MULTIPLICITY, THE NEUTRON CROSS
C       SECTION FOR MT=3 MUST BE ADJUSTED TO REPRESENT THE SAME
C       DATA AS MT=3 DOES IN THE PHOTON DATA
  240 ID=2
C       OBTAIN NEUTRON ELASTIC SCATTERING CROSS SECTION
      L1=LDICT(ID,IIN)
      IF(L1.EQ.0)GO TO 250
      LS1=IDICTS(ID,IIN)+LMOX2
      LEN=L1/2
      CALL TBSPLT(D(LS1),EOLD,LEN,XSIG2)
C       SUBTRACT THE ELASTIC SCATTERING CROSS SECTION FROM THE TOTAL
C       CROSS SECTION TO OBTAIN BASE NEUTRON MT=3 REACTION
      SIGMT3=SIGT-XSIG2
      GO TO 260
  250 SIGMT3=SIGT
  260 CONTINUE
C       SCAN THE PHOTON PRODUCTION DICTIONARY FOR ALL MT NUMBERS
C       NOT EQUAL TO MT=3
      DO 300 IX=1,NNTX
         MTG=LGCB(2*IX-1,IIN)
         IF(MTG.EQ.3)GO TO 300
         L1=LGCB2(2*IX,IIN)
         IF(L1.EQ.0)GO TO 370
         LS1=IGCBS2(2*IX,IIN)+LMOX4
         LEN=L1/2
         CALL TBSPLT(D(LS1),EOLD,LEN,SIGEX)
C     IF THE TOTAL PHOTON PRODUCTION CROSS SECTION IS ZERO AT
C     THE NEUTRON ENERGY, THEN THE NEUTRON CROSS SECTION SHOULD
C     NOT BE SUBTRACTED FROM MT3 TO MAINTAIN PROPER NORMALIZATION
         IF(SIGEX.EQ.0.0)GO TO 300
C     SET THE NEUTRON DICTIONARY ID NUMBER CORRESPONDING TO MTG
         IF((MTG.LT.51).OR.(MTG.GT.91))GO TO 270
         ID=14
         IMT3=MTG-51
         ID=ID+IMT3
  270    IF(MTG.EQ.4)ID=3
         IF(MTG.EQ.16)ID=8
         IF(MTG.EQ.17)ID=9
         IF(MTG.EQ.18)ID=10
         IF(MTG.EQ.22)ID=11
         IF(MTG.EQ.24)ID=12
         IF(MTG.EQ.28)ID=13
         IF(MTG.EQ.102)ID=55
         IF(MTG.EQ.103)ID=56
         IF(MTG.EQ.104)ID=57
         IF(MTG.EQ.105)ID=58
         IF(MTG.EQ.106)ID=59
         IF(MTG.EQ.107)ID=60
         IF(MTG.EQ.108)ID=61
         IF(MTG.EQ.109)ID=62
         IF(MTG.EQ.111)ID=63
         IF(MTG.EQ.112)ID=64
         IF(MTG.EQ.113)ID=65
         IF(MTG.EQ.114)ID=66
C     OBTAIN THE NEUTRON CROSS SECTION CORRESPONDING TO MTG AND
C     SUBTRACT IT OFF OF THE BASE NEUTRON MT=3 CROSS SECTION
         L1=LDICT(ID,IIN)
         IF(L1.EQ.0)GO TO 280
         LS1=IDICTS(ID,IIN)+LMOX2
         LEN=L1/2
         CALL TBSPLT(D(LS1),EOLD,LEN,XSIG)
         GO TO 290
  280    XSIG=0.0
  290    SIGMT3=SIGMT3-XSIG
         IF(SIGMT3.LE.0.0)GO TO 310
  300 CONTINUE
C    DETERMINE THE PHOTON MULTIPLICITY (YP)
      YP=SIGP/SIGMT3
      IF(YP.GE.100.0)GO TO 310
      GO TO 330
  310 CONTINUE
C       THIS SECTION OF CODING IS INCLUDED TO ACCOUNT FOR ANY
C       ENDF/B DATA INCONSISTENCY WHICH COULD YIELD A PHOTON OF
C       CONSIDERABLE WEIGHT.  THE FOLLOWING CODING WILL SAMPLE THE
C       PHOTON WEIGHT FROM THE GENERAL PHOTON YIELD ARRAY AND
C       ADJUST THE WEIGHT TO PHOTONS PER NON-ELASTIC COLLISION.
      L1=LGAM(IIN)
      IF(L1.EQ.0)GO TO 320
      LS1=IGAMS(IIN)+LMOX2
      LEN=L1/2
      CALL TBSPLT(D(LS1),EOLD,LEN,YP)
      YP=(YP*SIGT)/(SIGT-XSIG2)
      GO TO 330
  320 YP=1.00
C       THE FOLLOWING SECTION OF CODING IS INCLUDED TO DISTRIBUTE
C       THE WEIGHT ENDF/B-V DATA MAY GIVE A PARTICULAR PHOTON.
C       FOR EXAMPLE, ENDF/B-V DATA MAY ASSIGN A MULITPLICITY OF
C       75 TO A PARTICULAR PHOTON.  BECAUSE SUCH A PHOTON COULD
C       CONSIDERABLY MODIFY THE RESULTS OF A DETECTOR RESPONSE, THE
C       MULTIPLICITY (PHOTON WEIGHT) IS DISTRIBUTED TO SEVERAL
C       PHOTONS (SPLITTING OF SORTS) WITH BOTH WEIGHT AND ENERGY
C       BEING CONSERVED.  THIS RARELY OCCURS BUT IS NECESSARY.
  330 CONTINUE
C poisson distributed photon multiplicity CZ 13.8.92
      IGTRY=0
      MGPAR=INT(FLOAT(MAXPAR)*0.7)
  340 CALL G3POISS(YP,NUMBG,1)
      IGTRY=IGTRY+1
      IF(NUMBG.GT.INT(4.*YP).OR.
     +   NUMBG.GT.MGPAR.AND.IGTRY.LT.5) GOTO 340
      NUMBG=MIN(NUMBG,MGPAR)
C Allow 0 Photond to be generated
      IF(NUMBG.EQ.0) RETURN
      EPTOT = YP*EP
      EPSUM = 0.0
      DO 350 I=1,NUMBG
C       ASSUME ISOTROPIC PHOTON EMISSION IN THE LABORATORY SYSTEM
         CALL GTISO(U1,V1,W1)
C       SET THE PHOTON EXIT PARAMETERS
         UP=U1
         VP=V1
         WP=W1
         AGEP=AGE
         MTP=MT
C re-sample photon energy depending on model used CZ 13.8.92
         IF(IEP.EQ.2) THEN
            CALL PARTXS(D(LS2),D(LS2),EOLD,SIGP,EP1)
            IF(EP1.GT.0.0) EP=EP1
         ENDIF
         IF(IEP.EQ.1) THEN
            CALL SECEGY(EP1,D(LS3),EOLD,D(LS3))
            IF(EP1.GT.0.0) EP=EP1
         ENDIF
         EPSUM = EPSUM+EP
C check for energy conservation
         IF(EPSUM.GT.EPTOT.OR.I.EQ.NUMBG) EP = EPTOT-EPSUM+EP
C       STORE THE PHOTON
         CALL STOPAR(IDGAMA,NGAMA)
C end photon production when energy is used up  CZ 13.8.92
         IF(EPSUM.GT.EPTOT) GOTO 360
  350 CONTINUE
  360 RETURN
  370 WRITE(IOUT,10000)
10000 FORMAT(' PHOTON: THE PHOTON PRODUCTION ',
     +       'CROSS SECTION DATA WAS NOT FOUND (L1=0)')
      GOTO 390
  380 WRITE(IOUT,10100)
10100 FORMAT(' PHOTON: NO SECONDARY ENERGY ',
     +   'DISTRIBUTION WAS FOUND FOR THE CONTINUUM REACTION CHOSEN')
  390 WRITE(6,*) ' CALOR: ERROR in PHOTON ===> STOP '
      STOP
      END
